1 Define Parameters

# Paths
path_to_save_obj <- "../results"
path_to_saved_post_QC_obj <- str_c(
  path_to_save_obj,
  "seurat_object_cite_seq_seurat_doublets_excluded.rds",
  sep = "/"
)
path_to_save_seurat_wnn_obj <- str_c(
  path_to_save_obj,
  "seurat_object_cite_seq_seurat_wnn.rds",
  sep = "/"
)
saved_cell_cycle_obj <- str_c(
  "../data/cycle.rda",
  sep = "/"
)

1.1 Read saved seurat object post QC

filtered_bcll.combined <- readRDS(path_to_saved_post_QC_obj)

2 ADT and RNA multimodal analysis

DefaultAssay(filtered_bcll.combined) <- 'RNA'
filtered_bcll.combined <- NormalizeData(filtered_bcll.combined) %>% FindVariableFeatures() %>% ScaleData() %>% RunPCA()  %>%
  RunHarmony(reduction = "pca",dims = 1:30,group.by.vars = "gemid",assay.use = "RNA",reduction.save = "harmony_RNA")
## Centering and scaling data matrix
## PC_ 1 
## Positive:  MALAT1, IL32, CD3D, FYB1, GIMAP7, GIMAP4, GIMAP5, SRGN, CCR7, IL7R 
##     LCP2, CD3G, S100A4, TCF7, CD2, CD5, ITM2A, PLAAT4, SLFN5, TC2N 
##     CD69, INPP4B, ITK, TRBC1, SIRPG, GBP2, CD28, LAT, AHNAK, LINC00861 
## Negative:  STMN1, HMGB2, MKI67, NUSAP1, TOP2A, TUBA1B, HIST1H4C, HIST1H1B, H2AFZ, CENPF 
##     CDK1, HMGN2, HMGB1, TUBB, KIFC1, TUBB4B, KNL1, CCNB2, TPX2, BIRC5 
##     HIST1H3B, CCNB1, MCM7, CDCA8, HJURP, CDCA3, PTTG1, NCAPH, CENPM, AURKB 
## PC_ 2 
## Positive:  CD74, CD83, JCHAIN, FCRL3, FCRL5, IGKC, MZB1, LINC00926, DERL3, SPIB 
##     IGHG1, IGHM, FCRL2, IGHG3, FCGR2B, TNFRSF13B, CD1C, ENTPD1, SIGLEC10, FCRL4 
##     SERPINA9, CD24, IGLC2, IFNGR1, PLCG2, CD86, TNFRSF17, MYO1E, FAM30A, GRN 
## Negative:  IL32, CD3D, FYB1, GIMAP4, GIMAP7, GIMAP5, CD3G, CD2, LCP2, SRGN 
##     CD5, ITM2A, TCF7, IL7R, CD28, SLFN5, LAT, SIRPG, S100A4, TIGIT 
##     ITK, TC2N, CTLA4, INPP4B, FYN, TRBC1, CD247, GBP2, PLAAT4, CCR7 
## PC_ 3 
## Positive:  PCNA, MCM7, MCM3, GINS2, RRM2, CLSPN, MCM6, CDC6, FEN1, MCM4 
##     DHFR, HELLS, FAM111B, ATAD2, TK1, CDC45, DTL, RFC2, MCM10, CCNE2 
##     WDR76, UHRF1, FABP5, PAICS, RFC4, EXO1, PRKDC, LIG1, RRM1, HIST1H3C 
## Negative:  KIF20A, CENPE, CCNB2, HMMR, CENPF, AURKA, CCNB1, ASPM, KPNA2, ARL6IP1 
##     DEPDC1, DLGAP5, KIF14, PLK1, CDC20, TOP2A, NEK2, PSRC1, KNSTRN, KIF23 
##     CDCA8, SGO2, TPX2, UBE2C, BUB1, CKAP2, HJURP, NUF2, CDKN3, ECT2 
## PC_ 4 
## Positive:  LYZ, TYROBP, FCER1G, LGALS2, LGALS1, HSP90B1, FGL2, CPVL, HSPA5, NME2 
##     SSR4, JAML, CHCHD2, DERL3, GRN, ANXA2, IGSF6, S100A6, MZB1, GSN 
##     RPL8, SLAMF7, SEC11C, MS4A6A, SERPINA1, SUB1, CALR, S100A11, CRIP1, JCHAIN 
## Negative:  MTRNR2L12, MALAT1, HIST1H3C, HIST1H2BF, HIST1H3B, MTRNR2L8, HIST1H4F, HIST1H1B, HIST1H1D, HIST1H4D 
##     HIST1H2AH, HIST1H2BJ, HIST1H2BH, HIST1H1C, HIST1H3F, HIST1H2BB, HIST1H2BL, HIST2H2BF, HIST1H2AE, HIST1H4L 
##     HIST1H2BM, HIST1H2BC, HIST1H1E, HIST1H4H, HIST1H3H, HIST1H2BN, HIST1H2BI, HIST1H4B, HIST1H2AB, HIST1H2BO 
## PC_ 5 
## Positive:  JCHAIN, MZB1, DERL3, IGHG1, SSR4, XBP1, SEC11C, LINC02362, TENT5C, HSP90B1 
##     IGKC, IGHG4, FKBP11, IGHG3, AC012236.1, IGLC2, IGHGP, TXNDC5, SLAMF7, PRDX4 
##     LMAN1, IGLC1, TNFRSF17, IGHG2, HYOU1, IGLC3, PDIA4, RASSF6, HSPA5, ZBP1 
## Negative:  LYZ, ACTB, TYROBP, CAPG, FGL2, IFNGR1, LGALS2, FCER1G, CPVL, FCRL4 
##     MT-CO1, CD74, ZEB2, LST1, GSN, TNFRSF13B, IGSF6, CD83, CD1C, AIF1 
##     CCR6, ITGAX, PLAC8, MNDA, JAML, SERPINA1, HSP90AB1, SPIB, FGR, CD68
## Harmony 1/10
## Harmony 2/10
## Harmony 3/10
## Harmony 4/10
## Harmony 5/10
## Harmony converged after 5 iterations
## Warning: Invalid name supplied, making object name syntactically valid. New
## object name is Seurat..ProjectDim.RNA.harmony_RNA; see ?make.names for more
## details on syntax validity

3 Check for Cell cycle stages variability between clusters

load(saved_cell_cycle_obj)
filtered_bcll.combined <- CellCycleScoring(filtered_bcll.combined, 
                                 g2m.features = g2m_genes, 
                                 s.features = s_genes)
DefaultAssay(filtered_bcll.combined) <- 'ADT'
# we will use all ADT features for dimensional reduction
# we set a dimensional reduction name to avoid overwriting the 
VariableFeatures(filtered_bcll.combined) <- rownames(filtered_bcll.combined[["ADT"]])
filtered_bcll.combined <- NormalizeData(filtered_bcll.combined, normalization.method = 'CLR', margin = 2) %>% 
  ScaleData() %>% RunPCA(reduction.name = 'apca') %>%
  RunHarmony(reduction = "apca",dims = 1:20, group.by.vars = "gemid", assay.use = "ADT", reduction.save = "harmony_protein")
## Normalizing across cells
## Centering and scaling data matrix
## PC_ 1 
## Positive:  HLA-DR, CD138-(Syndecan-1), CD19.1, CD20, CD45RA, CD21, CD90-(Thy1), NLRP2.1, CD40.1, 0856-anti-Tau-Phospho-(Thr181) 
##     IgA, IgG-Fc, 1055-anti-c-Met, TCR-gamma-delta, HLA-F.1, CD38.1, CD268-(BAFF-R), CD133, CD81-(TAPA-1), HLA-A-B-C 
##     CD193-(CCR3), CD309-(VEGFR2), CD30, CD35, CD15-(SSEA-1), CD10, CD54, CD22.1, CD11b, CD47.1 
## Negative:  CD2.1, CD5.1, TCR-alpha-beta, CD4.1, CD3, CD49f, CD278-(ICOS), CD7.1, CD62L, CD45RO 
##     CD28.1, CD274-(B7-H1--PD-L1), CD66b, CD197-(CCR7), TIGIT-(VSTM3), CD223-(LAG-3), CD122-(IL-2Rbeta), CD127-(IL-7Ralpha), IgG2a--kappa-isotype-Ctrl, CD23 
##     CD224, CD279-(PD-1), CD267-(TACI), CD96-(TACTILE), CD57-Recombinant, CD303-(BDCA-2), CD169-(Sialoadhesin--Siglec-1), CD370-(CLEC9A-DNGR1), CD101-(BB27), CD106 
## PC_ 2 
## Positive:  Podoplanin, CD35, CX3CR1.1, CD70.1, CD58-(LFA-3), TCR-Vbeta13.1, CD336-(NKp44), TCR-Vgamma9, CD269-(BCMA), Mac-2-(Galectin-3) 
##     CD32, CD107a-(LAMP-1), TCR-Valpha24-Jalpha18-(iNKT-cell), CD21, CD82.1, CD326-(Ep-CAM), CD307d-(FcRL4), CD54, CD183-(CXCR3), CD360-(IL-21R) 
##     CD150-(SLAM), CD268-(BAFF-R), CD1c, CD73-(Ecto-5-nucleotidase), CD196-(CCR6), CD24.1, CD207.1, CD144-(VE-Cadherin), CD39, CD40.1 
## Negative:  CD2.1, CD3, CD7.1, CD5.1, TCR-alpha-beta, CD4.1, CD278-(ICOS), CD27.1, CD28.1, CD8 
##     CD11a, CD49f, CD279-(PD-1), CD18, CD45RO, CD45, IgA, IgG2b--kappa-isotype-Ctrl, CD99.1, CD62L 
##     CD127-(IL-7Ralpha), CD95-(Fas), HLA-A2, CD194-(CCR4), CD161, CD328-(Siglec-7), CD138-(Syndecan-1), CD36.1, CD38.1, CD90-(Thy1) 
## PC_ 3 
## Positive:  CD45, CD99.1, CD95-(Fas), CD52.1, CD11a, CD278-(ICOS), CD5.1, CD2.1, CD18, CD28.1 
##     CD47.1, CD4.1, CD3, CD69.1, CD45RO, CD27.1, CD82.1, CD279-(PD-1), TCR-alpha-beta, CD7.1 
##     CD38.1, CD224, CD44.1, CD62L, HLA-A-B-C, CD81-(TAPA-1), CD54, CD275-(B7-H2--ICOSL), CD39, CD194-(CCR4) 
## Negative:  CD303-(BDCA-2), CD169-(Sialoadhesin--Siglec-1), CD23, CD154, CD158f-(KIR2DL5), CD124-(IL-4Ralpha), CD204, CD66b, CD370-(CLEC9A-DNGR1), CD206-(MMR) 
##     CD269-(BCMA), CD163.1, CD122-(IL-2Rbeta), CD144-(VE-Cadherin), CD106, CD178-(Fas-L), CD209-(DC-SIGN), CD267-(TACI), CD141-(Thrombomodulin), CD34.1 
##     integrin-beta7, TSLPR-(TSLP-R), CD274-(B7-H1--PD-L1), TCR-Valpha24-Jalpha18-(iNKT-cell), Mac-2-(Galectin-3), CD294-(CRTH2), CD324-(E-Cadherin), IgG1--kappa-isotype-Ctrl, CD197-(CCR7), FcepsilonRIalpha 
## PC_ 4 
## Positive:  CD90-(Thy1), CD138-(Syndecan-1), 0856-anti-Tau-Phospho-(Thr181), 1055-anti-c-Met, NLRP2.1, CD226-(DNAM-1), HLA-F.1, CD133, B7-H4, CD154 
##     TCR-gamma-delta, IgG-Fc, CD309-(VEGFR2), CD15-(SSEA-1), CD193-(CCR3), CD279-(PD-1), CD337-(NKp30), CD30, CD336-(NKp44), CD45RO 
##     CD137L-(4-1BB-Ligand), IgA, CD58-(LFA-3), CD158f-(KIR2DL5), CD269-(BCMA), XCR1.1, CD28.1, CD303-(BDCA-2), CD124-(IL-4Ralpha), CD278-(ICOS) 
## Negative:  CD32, CD73-(Ecto-5-nucleotidase), CD35, CD1c, CD196-(CCR6), IgD, CD39, CD305-(LAIR1), CD21, IgM 
##     CD268-(BAFF-R), CD82.1, CD54, CD22.1, CD45RA, CD19.1, CD31, CD49d, CD71, CD47.1 
##     CD44.1, CD23, CD69.1, Ig-light-chain-lambda, CD52.1, CD62L, CD123, IgG2a--kappa-isotype-Ctrl, CD40.1, CD307d-(FcRL4) 
## PC_ 5 
## Positive:  IgD, CD305-(LAIR1), IgM, CD69.1, CD44.1, CD196-(CCR6), CD49d, CD1c, CD90-(Thy1), CD24.1 
##     CD73-(Ecto-5-nucleotidase), CD138-(Syndecan-1), CD32, 0856-anti-Tau-Phospho-(Thr181), CD11b, CD31, 1055-anti-c-Met, TCR-gamma-delta, NLRP2.1, CD133 
##     CD35, HLA-F.1, CD15-(SSEA-1), CD309-(VEGFR2), CD30, CD39, IgG-Fc, Ig-light-chain-lambda, CD193-(CCR3), B7-H4 
## Negative:  CD10, CD81-(TAPA-1), CD71, CD38.1, CD20, CD95-(Fas), CD54, CD146, CD40.1, HLA-A-B-C 
##     CD19.1, CD86.1, CD98, CD80.1, CD45, CD58-(LFA-3), HLA-DR, IgG2a--kappa-isotype-Ctrl, CD184-(CXCR4), CD57-Recombinant 
##     CD197-(CCR7), CD124-(IL-4Ralpha), CD274-(B7-H1--PD-L1), CD107a-(LAMP-1), CD82.1, CD319-(CRACC), CD66b, CD106, CD370-(CLEC9A-DNGR1), TIGIT-(VSTM3)
## Warning: Cannot add objects with duplicate keys (offending key: PC_), setting
## key to 'apca_'
## Harmony 1/10
## Harmony 2/10
## Harmony converged after 2 iterations
## Warning: Invalid name supplied, making object name syntactically valid. New
## object name is Seurat..ProjectDim.ADT.harmony_protein; see ?make.names for more
## details on syntax validity

4 Multimodal Neighbor Identification

# Identify multimodal neighbors. These will be stored in the neighbors slot, 
# and can be accessed using filtered_bcll.combined[['weighted.nn']]
# The WNN graph can be accessed at filtered_bcll.combined[["wknn"]], 
# and the SNN graph used for clustering at filtered_bcll.combined[["wsnn"]]
# Cell-specific modality weights can be accessed at filtered_bcll.combined$RNA.weight
filtered_bcll.combined <- FindMultiModalNeighbors(
  filtered_bcll.combined, reduction.list = list("harmony_RNA", "harmony_protein"),
  dims.list = list(1:30, 1:20), modality.weight.name = "RNA.weight"
)
## Calculating cell-specific modality weights
## Finding 20 nearest neighbors for each modality.
## Calculating kernel bandwidths
## Warning in FindMultiModalNeighbors(filtered_bcll.combined, reduction.list
## = list("harmony_RNA", : The number of provided modality.weight.name is not
## equal to the number of modalities. RNA.weight ADT.weight are used to store the
## modality weights
## Finding multimodal neighbors
## Constructing multimodal KNN graph
## Constructing multimodal SNN graph

5 UMAP visualisation

feat_clust <- c(0.1,0.5,1,1.5,2)
filtered_bcll.combined <- RunUMAP(filtered_bcll.combined, nn.name = "weighted.nn", reduction.name = "wnn.umap", reduction.key = "wnnUMAP_")
## Warning: The default method for RunUMAP has changed from calling Python UMAP via reticulate to the R-native UWOT using the cosine metric
## To use Python UMAP via reticulate, set umap.method to 'umap-learn' and metric to 'correlation'
## This message will be shown once per session
## 17:07:18 UMAP embedding parameters a = 0.9922 b = 1.112
## 17:07:19 Commencing smooth kNN distance calibration using 1 thread
## 17:07:20 Initializing from normalized Laplacian + noise
## 17:07:20 Commencing optimization for 200 epochs, with 1424488 positive edges
## 17:07:37 Optimization finished
filtered_bcll.combined <- FindClusters(filtered_bcll.combined, graph.name = "wsnn", algorithm = 3, resolution = feat_clust, verbose = FALSE)
Seurat::DimPlot(filtered_bcll.combined, reduction = 'wnn.umap', label = TRUE, repel = TRUE, label.size = 2.5, group.by = c("wsnn_res.0.1","wsnn_res.0.5","wsnn_res.1","wsnn_res.1.5","wsnn_res.2")) + NoLegend()  & Seurat::NoLegend()
## Warning: ggrepel: 2 unlabeled data points (too many overlaps). Consider
## increasing max.overlaps
## Warning: ggrepel: 3 unlabeled data points (too many overlaps). Consider
## increasing max.overlaps

filtered_bcll.combined <- FindClusters(filtered_bcll.combined, graph.name = "wsnn", algorithm = 3, resolution = 0.1, verbose = FALSE)
filtered_bcll.combined <- RunUMAP(filtered_bcll.combined, reduction = 'pca', dims = 1:30, assay = 'RNA', 
              reduction.name = 'rna.umap', reduction.key = 'rnaUMAP_')
## 17:11:14 UMAP embedding parameters a = 0.9922 b = 1.112
## 17:11:14 Read 42929 rows and found 30 numeric columns
## 17:11:14 Using Annoy for neighbor search, n_neighbors = 30
## 17:11:14 Building Annoy index with metric = cosine, n_trees = 50
## 0%   10   20   30   40   50   60   70   80   90   100%
## [----|----|----|----|----|----|----|----|----|----|
## **************************************************|
## 17:11:18 Writing NN index file to temp file /tmp/Rtmpordw26/file105cb4571c3d33
## 17:11:18 Searching Annoy index using 1 thread, search_k = 3000
## 17:11:31 Annoy recall = 100%
## 17:11:31 Commencing smooth kNN distance calibration using 1 thread
## 17:11:33 Initializing from normalized Laplacian + noise
## 17:11:34 Commencing optimization for 200 epochs, with 1935756 positive edges
## 17:11:49 Optimization finished
filtered_bcll.combined <- RunUMAP(filtered_bcll.combined, reduction = 'apca', dims = 1:18, assay = 'ADT', 
              reduction.name = 'adt.umap', reduction.key = 'adtUMAP_')
## 17:11:49 UMAP embedding parameters a = 0.9922 b = 1.112
## 17:11:49 Read 42929 rows and found 18 numeric columns
## 17:11:49 Using Annoy for neighbor search, n_neighbors = 30
## 17:11:49 Building Annoy index with metric = cosine, n_trees = 50
## 0%   10   20   30   40   50   60   70   80   90   100%
## [----|----|----|----|----|----|----|----|----|----|
## **************************************************|
## 17:11:53 Writing NN index file to temp file /tmp/Rtmpordw26/file105cb47044f1e6
## 17:11:53 Searching Annoy index using 1 thread, search_k = 3000
## 17:12:05 Annoy recall = 100%
## 17:12:05 Commencing smooth kNN distance calibration using 1 thread
## 17:12:06 Initializing from normalized Laplacian + noise
## 17:12:07 Commencing optimization for 200 epochs, with 1871324 positive edges
## 17:12:22 Optimization finished
p1 <- DimPlot(filtered_bcll.combined, reduction = 'rna.umap', label = TRUE, 
              repel = TRUE, label.size = 2.5) + NoLegend()
p2 <- DimPlot(filtered_bcll.combined, reduction = 'adt.umap', label = TRUE, 
              repel = TRUE, label.size = 2.5) + NoLegend()
p1

p2

Seurat::DimPlot(filtered_bcll.combined, reduction = 'wnn.umap', label = TRUE, repel = TRUE, label.size = 2.5, group.by = "wsnn_res.0.1")

metadata <- filtered_bcll.combined@meta.data
df = count(metadata,metadata$wsnn_res.0.1)
colnames(df) = c("Cluster","Number_of_cells")
library(knitr)
library(kableExtra)
## 
## Attaching package: 'kableExtra'
## The following object is masked from 'package:dplyr':
## 
##     group_rows
kable(df) %>%
  kable_styling("striped", full_width = T)
Cluster Number_of_cells
0 10590
1 10446
10 158
11 2
2 9562
3 5851
4 2156
5 1620
6 989
7 735
8 461
9 359

6 Check Batch effect

6.1 Check for subproject variability

DimPlot(filtered_bcll.combined, reduction = "wnn.umap", group.by = 'subproject')

6.2 Cell Cycle

# Plot the PCA colored by cell cycle phase
DimPlot(filtered_bcll.combined,
        reduction = "pca",
        group.by= "Phase",
        split.by = "Phase")

DimPlot(filtered_bcll.combined, reduction = "wnn.umap", group.by = 'Phase')

6.3 Scrublet Doublets

DimPlot(filtered_bcll.combined, reduction = "wnn.umap", group.by = 'scrublet_predicted_doublet')

FeaturePlot(filtered_bcll.combined, features = "scrublet_doublet_scores")

6.4 Check vireo doublets

DimPlot(filtered_bcll.combined, group.by = "genotype_based_doublet_flag", pt.size = 0.1)

6.5 Unassigned to any donor

DimPlot(filtered_bcll.combined, group.by = "genotype_based_unassigned_flag", pt.size = 0.1)

6.6 Plot UMAP with respect to TCR

DimPlot(filtered_bcll.combined, group.by = "tcr_flag", pt.size = 0.1)

6.7 Plot UMAP with respect to BCR

DimPlot(filtered_bcll.combined, group.by = "bcr_flag", pt.size = 0.1)

6.8 Plot MAIT cells

DimPlot(filtered_bcll.combined, group.by = "mait_evidence", pt.size = 0.1)

6.9 Plot iNKT cells

DimPlot(filtered_bcll.combined, group.by = "inkt_evidence", pt.size = 0.1)

7 RNA and ADT weight per cluster

7.1 Violin plot to check the RNA weight

Seurat::VlnPlot(filtered_bcll.combined, features = "RNA.weight", sort = TRUE, pt.size = 0.1) & Seurat::NoLegend()

7.2 Violin plot to check the ADT weight

Seurat::VlnPlot(filtered_bcll.combined, features = "ADT.weight", sort = TRUE, pt.size = 0.1) & Seurat::NoLegend()

7.3 Save the seurat object

saveRDS(filtered_bcll.combined, file = path_to_save_seurat_wnn_obj)

8 Session Info

sessionInfo()
## R version 4.0.3 (2020-10-10)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 20.04.2 LTS
## 
## Matrix products: default
## BLAS:   /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.9.0
## LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.9.0
## 
## locale:
##  [1] LC_CTYPE=en_GB.UTF-8       LC_NUMERIC=C              
##  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_GB.UTF-8    
##  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_GB.UTF-8   
##  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
##  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] kableExtra_1.3.4   harmony_0.1.0      Rcpp_1.0.8.3       forcats_0.5.1     
##  [5] stringr_1.4.0      purrr_0.3.4        readr_2.1.2        tidyr_1.2.0       
##  [9] tibble_3.1.6       tidyverse_1.3.1    scales_1.1.1       patchwork_1.1.1   
## [13] ggplot2_3.3.5      dplyr_1.0.8        SeuratObject_4.0.4 Seurat_4.1.0      
## [17] knitr_1.37        
## 
## loaded via a namespace (and not attached):
##   [1] readxl_1.3.1          backports_1.4.1       systemfonts_1.0.4    
##   [4] plyr_1.8.6            igraph_1.2.11         lazyeval_0.2.2       
##   [7] splines_4.0.3         listenv_0.8.0         scattermore_0.8      
##  [10] digest_0.6.29         htmltools_0.5.2       fansi_1.0.2          
##  [13] magrittr_2.0.2        tensor_1.5            cluster_2.1.2        
##  [16] ROCR_1.0-11           tzdb_0.2.0            globals_0.14.0       
##  [19] modelr_0.1.8          matrixStats_0.61.0    svglite_2.1.0        
##  [22] spatstat.sparse_2.1-0 colorspace_2.0-3      rvest_1.0.2          
##  [25] ggrepel_0.9.1         haven_2.4.3           xfun_0.30            
##  [28] crayon_1.5.0          jsonlite_1.8.0        spatstat.data_2.1-2  
##  [31] survival_3.3-1        zoo_1.8-9             glue_1.6.2           
##  [34] polyclip_1.10-0       gtable_0.3.0          webshot_0.5.2        
##  [37] leiden_0.3.9          future.apply_1.8.1    abind_1.4-5          
##  [40] DBI_1.1.2             spatstat.random_2.1-0 miniUI_0.1.1.1       
##  [43] viridisLite_0.4.0     xtable_1.8-4          reticulate_1.24      
##  [46] spatstat.core_2.4-0   htmlwidgets_1.5.4     httr_1.4.2           
##  [49] RColorBrewer_1.1-2    ellipsis_0.3.2        ica_1.0-2            
##  [52] farver_2.1.0          pkgconfig_2.0.3       sass_0.4.0           
##  [55] uwot_0.1.11           dbplyr_2.1.1          deldir_1.0-6         
##  [58] utf8_1.2.2            tidyselect_1.1.2      labeling_0.4.2       
##  [61] rlang_1.0.2           reshape2_1.4.4        later_1.3.0          
##  [64] munsell_0.5.0         cellranger_1.1.0      tools_4.0.3          
##  [67] cli_3.2.0             generics_0.1.2        broom_0.7.12         
##  [70] ggridges_0.5.3        evaluate_0.15         fastmap_1.1.0        
##  [73] yaml_2.3.5            goftest_1.2-3         fs_1.5.2             
##  [76] fitdistrplus_1.1-8    RANN_2.6.1            pbapply_1.5-0        
##  [79] future_1.24.0         nlme_3.1-155          mime_0.12            
##  [82] xml2_1.3.3            compiler_4.0.3        rstudioapi_0.13      
##  [85] plotly_4.10.0         png_0.1-7             spatstat.utils_2.3-0 
##  [88] reprex_2.0.1          bslib_0.3.1           stringi_1.7.6        
##  [91] highr_0.9             RSpectra_0.16-0       lattice_0.20-45      
##  [94] Matrix_1.4-0          vctrs_0.3.8           pillar_1.7.0         
##  [97] lifecycle_1.0.1       spatstat.geom_2.3-2   lmtest_0.9-39        
## [100] jquerylib_0.1.4       RcppAnnoy_0.0.19      data.table_1.14.2    
## [103] cowplot_1.1.1         irlba_2.3.5           httpuv_1.6.5         
## [106] R6_2.5.1              promises_1.2.0.1      KernSmooth_2.23-20   
## [109] gridExtra_2.3         parallelly_1.30.0     codetools_0.2-18     
## [112] MASS_7.3-55           assertthat_0.2.1      withr_2.5.0          
## [115] sctransform_0.3.3     mgcv_1.8-39           parallel_4.0.3       
## [118] hms_1.1.1             grid_4.0.3            rpart_4.1.16         
## [121] rmarkdown_2.13        Rtsne_0.15            shiny_1.7.1          
## [124] lubridate_1.8.0